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It has previously been shown that the network of connected minima on a potential energy land- 
scape is scale-free, and that this reflects a power-law distribution for the areas of the basins of 
attraction surrounding the minima. Here, we set out to understand more about the physical origins 
of these puzzling properties by examining how the potential energy landscape of a 13-atom cluster 
evolves with the range of the potential. In particular, on decreasing the range of the potential the 
number of stationary points increases and thus the landscape becomes rougher and the network gets 
larger. Thus, we are able to follow the evolution of the potential energy landscape from one with 
just a single minimum to a complex landscape with many minima and a scale-free pattern of con- 
nections. We find that during this growth process, new edges in the network of connected minima 
preferentially attach to more highly-connected minima, thus leading to the scale-free character. Fur- 
thermore, minima that appear when the range of the potential is shorter and the network is larger 
have smaller basins of attraction. As there are many of these smaller basins because the network 
grows exponentially, the observed growth process thus also gives rise to a power-law distribution for 
the hyperareas of the basins. 



I. INTRODUCTION 

A potential energy landscape (PEL) is a high- 
dimensional surface describing how the potential energy 
of a system varies with the coordinates of each atom in 
the system^ Since PELs control the structure, thermo- 
dynamics and dynamics of a system, there has been an 
ongoing research programme that has sought to under- 
stand how the behaviour of a system can be related back 
to features of a complex PEL. This approach has led to 
new insights into protein folding^ and the nature of the 
glass transition^ In these applications the emphasis is 
on how PELs differ, for example between proteins that 
are good and bad folders or liquids that are fragile and 
strong. In this paper by contrast, we are interested in 
more universal properties of high dimensional PELs. 

When characterizing PELs the focus is often on the 
stationary points of these surfaces. Some fundamental 
properties associated with stationary points of PELs are 
well established. For example, the number of station- 
ary points increases exponentially with system size4^£ 
The energy distribution of both minima and transition 
states follows a Gaussian distribution, as has been seen 
empirically^ and predicted theoretically from the cen- 
tral limit theorem. 9 More recently, network-derived in- 
sights have been obtained for PELsr^^^^^ To ob- 
tain a network representation of a PEL, each minima is 
mapped to a node in the network and each transition 
state to an edge, as in Fig. [1] The resulting "inher- 
ent structure" networ k 10 ^ 1 provides a dynamically moti- 
vated description of the connectivity and organization of 
the PEL, since at sufficiently low temperatures, a system 
vibrates within the basin of attraction of a minimum, 
or inherent structure^ before moving to another via a 




FIG. 1: (Colour online), (a) A model two-dimensional energy 
surface, (b) A contour plot of this surface illustrating the 
inherent structure partition of the configuration space into 
basins of attraction surrounding minima. The basin bound- 
aries are represented by thick lines, and the minima and tran- 
sition states by dots, (c) The resulting representation of the 
landscape as a network. 



transition state valleyJ^ 

For small Lennard- Jones clusters, these networks have 
been shown to have both scale-free and small world 



properties 
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The characteristic of a scale-free network 



is a power-law tail to the degree distribution^ where 
the degree of a node is the number of edges connected 
to it, and the characteristic of a small- world network is 
a logarithmic increase in the average separation between 
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nodes with network sizeJ^ These two properties are re- 
lated since the high-degree hubs present in scale-free net- 
works connect up the network, making the path lengths 
between nodes very short ^ Recently, a wide variety of 
networks have been shown to be scale free ; 20 ' 21 including 
the world wide web , 22 ' 23 social network o 24 ' 25 and biolog- 
ical networks i 26 i 27 Typically, the exponent of the power- 
law lies between —2 and —3, as is the case for the inherent 
structure networks, which have an exponent of approx- 
imately — 2.78J 10 ' 11 Interestingly, networks based on the 
configuration space of a protein ; 28 ' 29 although defined in 
a more coarse-grained way than the inherent structure 
networks, are also scale-free, suggesting that this feature 
could be a universal property of PELs. 

A popular explanation for the scale-free character of 
many of these networks is in terms of preferential at- 
tachment during the growth of these networks.— In pref- 
erential attachment, new nodes added to the network 
are more likely to attach to old nodes with high degree. 
Therefore, the degree of high-degree nodes increases at a 
faster rate, a 'rich get richer' growth process, and leads 
to a power-law degree distribution. 

By contrast, a PEL is static — it is just determined by 
the form of the potential and the system size. Instead, 
for PELs, the scale-free nature seems to be related to 
the distribution of basin areas. Each minimum on the 
PEL is surrounded by a basin of attraction, a region of 
configuration space where following the gradient downhill 
will lead to the same minimum. Low-energy minima have 
large basins of attraction ; 30 ' 31 and so unsurprisingly they 
have high degree because they can fit many neighbours 
along their long boundaries. A broad degree distribution 
thus reflects a broad basin area distribution. 

This connection between area and degree is more ex- 
plicit in Apollonian networks, which are model spa- 
tial, scale-free network o 32 ' 33 that have similar topologi- 
cal properties to inherent structure networks i 14 ' 33 These 
Apollonian networks are based on Apollonian packings 
in which all of space is tiled by hyperspheres (see Fig. [2] 
for a two-dimensional example) . Each hypersphere gives 
rise to a node in the network, and two nodes are con- 
nected if the corresponding hyperspheres touch. The 
packing is constructed iteratively and at each step in the 
generation of the packing each interstice is filled by the 
largest possible hypersphere that just touches the sur- 
rounding hyperspheres. The fractal, self-similar charac- 
ter of the packings leads both to the scale-free topology 
of the Apollonian networks, and a power-law distribution 
for the hypervolumes of the hyperspheres. 

The Apollonian packings and networks suggest how 
the basins of attraction on a PEL might have to tile con- 
figuration space in order to generate a scale-free inher- 
ent structure network. Namely, the tiling would need to 
be hierarchical with larger basins surrounded by smaller 
basins, which are in turn surrounded by smaller basins, 
and so on. Interestingly, when the basin area distribution 
has been probed for model liquids, a power-law distribu- 
tion has been found with the exponent expected from 




FIG. 2: A two-dimensional Apollonian packing. Space is 
filled with different sized disks, starting from an initial con- 
figuration where the three larger disks are placed within the 
bounding circle. The packing is generated iteratively, with 
the largest disk possible added to each gap at each iteration. 
This process is continued ad infinitum, thus filling the space 
with successively smaller and smaller disks. 



analogy to the Apollonian packings,— thus suggesting 
that the inherent structure networks for these systems 
are scale-free (note that it is effectively impossible to 
obtain these networks directly due to the high dimen- 
sionality of their configuration space). Although these 
results and ideas show how the topology of the inherent 
structure networks can be understood in terms of how 
the basins of attraction tile configuration space, why this 
tiling should be hierarchical, fractal and Apollonian-like 
is still a puzzle. 

Recently, using a simple model we explored whether 
these features could simply be driven by the variation in 
the energies of the minimal^ 4 - In this 'eggbox' model, we 
took the Gaussian distribution of energies of the minima 
as a starting point, and examined how the PEL landscape 
evolved (e.g. through the swallowing up of higher energy 
minima by deep basins) as the variance of the Gaussian 
increased. Although the resulting area distribution was 
considerably broader than for a flat landscape where all 
the minima have the same energy, these changes were not 
sufficient to genearate a power law tail for the degree or 
basin area distributions. 

Here, we provide an alternative perspective on why the 
inherent structure networks for PELs are scale-free, and, 
as we will see, hence obtain insights into the nature of the 
division of configuration space into basins of attraction. 
Although, as already mentioned, a PEL is static, there 
may be ways to envisage the generation of the network as 
a quasi-growth process, and hence to apply the ideas of 
preferential attachment. Recently, we looked at the role 
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FIG. 3: The generalized Lennard- Jones potential (Eq.[TJ for 
p = 2, 4 and 6. . 



played by the energetic ordering of minima. 13 Here, we 
wish to look at how the PEL (and the inherent structure 
network) evolves on going from a simple landscape with 
one minimum to a complex one with many (JV) minima 
and a scale-free topology. The reverse of this process, 
smoothing a complex landscape so that there are signif- 
icantly fewer minima, and perhaps only one, has been 
a common approach in global optimization . 35 ' 36 i 37 ' 38 i 39 
The idea behind this approach is that global optimiza- 
tion is likely to be much easier on the smoothed land- 
scape, and trivially so if there is only one minimum left. 
Of this set of global optimization algorithms, particularly 
relevant to us are those that achieve this smoothing by 
changing a potential parameter i 36 i 37 ' 39 

Here, we examine the evolution of the landscape of the 
13-atom cluster for a generalized Lennard-Jones potential 
with parameter p, defined by 



V = 4e 



r J V r 



(1) 



where r is the distance between two atoms, e is the well- 
depth and 2 1 / p o is the equilibrium pair separation. This 
is summed over all pairs of atoms to give the total poten- 
tial energy of a configuration. The range of the potential 
is determined by the parameter p, where for p = 6, the 
potential is equivalent to the Lennard-Jones potential. 
As p decreases, the repulsion becomes less steep and the 
attraction has a longer range. The net effect is that the 
width of the well increases, as shown in Fig. [3] 

The effects of p (or the analogous parameter for 
the Morse potential) on the landscape have been well 
stvdied&dMLg&M&SMUZ At sufficiently small p, 
there is a single minimum on the landscape. As p in- 
creases, the number of minima increases rapidly, because 
potential wells surrounding the minima become narrower, 
allowing new minima to appear at higher energy. Equiv- 
alently, as p decreases, the lower energy wells become 



wider, thus swallowing up the higher energy minima. 
Downhill barriers also increase as p increases. The net ef- 
fect is that a smooth landscape at low p becomes rougher. 

In Section HH we provide details of how the appear- 
ance and disappearance of minima and transition states 
can be described by catastrophe theor y 49 ! 50 and the ef- 
fect on the inherent structure network. In Section IITT1 we 
describe the method used to detect preferential attach- 
ment. In Section [TVl we discuss how the landscape for 
the 13-atom Lennard-Jones cluster evolves, and whether 
preferential attachment provides a useful description of 
that evolution. 



II. CATASTROPHES 

The two most common types of catastrophe for PELs 
are the fold and the cusp, illustrated in Fig. 2] In the fold, 
a minimum and a maximum move closer together as some 
parameter is varied, in this case as the range is increased. 
At the catastrophe the two stationary points meet and 
disappear, leaving just a shoulder in the landscape. The 
cusp in Fig. Hfb) involves two symmetry-related maxima 
and one minimum. Upon increasing the range, they move 
closer together until they collide, leaving one maximum. 

Fig. 3] is a one-dimensional picture, whereas the PEL 
is high-dimensional. However, these catastrophes can 
provide a good description of the evolution of one- 
dimensional 'reaction coordinates' that link stationary 
points on the PEL. 4 - These reaction coordinates are de- 
fined in terms of steepest-descent pathways, i.e. the paths 
follow the gradient downhill at every point, and thus 
there is no component of the gradient in the space orthog- 
onal to the coordinate. Hence the points where there is 
a maximum or minimum in the reaction coordinate cor- 
respond to stationary points of the PEL. 

We are interested in catastrophes involving minima 
and transition states, where minima are stationary points 
that are minimal in all directions, and transition states 
are maximal in one direction and minimal in all others, 
i.e. the Hessian has one negative eigenvalue. This leaves 
two possible scenarios. In the first, the minima in Fig.0] 
represent minima and the maxima transition states. In 
the second, the minima in Fig. U represent transition 
states and the maxima index 2 saddles, which are max- 
imal in one more dimension than transition states. The 
cusp catastrophe can only occur when there is symmetry 
present, as the pathways either side of the central sta- 
tionary point are equivalent. Therefore, any stationary 
points connected by the pathways must be permutational 
isomers, where the system has the same structure but 
some atoms swap positions. The low p maximum must 
therefore connect two permutational isomers of the same 
minimum. 

There are some subtleties in the way the inherent struc- 
ture network has been defined. Edges correspond to tran- 
sition states, but not all transition states give rise to 
edges. More than one transition state may connect the 
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(a) 



(b) 




FIG. 4: (Colour online) Changes in a landscape close to a 
catastrophe. The loci of the stationary points are also plot- 
ted, (a) The fold catastrophe. Initially, there are two sta- 
tionary points, a maximum and a minimum. Upon varying a 
parameter, decreasing p (increasing the range) in the case of 
the PEL, these stationary points move closer together until 
they collide and disappear, (b) The cusp catastrophe. Ini- 
tially, there is one minimum and two maxima, which collide 
upon increasing the range to leave a single maximum. 



same pair of nodes, forming a multiple edge (ME), and 
some connect permutational isomers of the same mini- 
mum, forming a self-connection (SC). ME/SCs are not 
considered in the inherent structure network as we are 
interested only in whether two minima are connected, 
not how many times, and we do not distinguish between 
permutational isomers of a minimum. This choice is im- 
portant for the network topologyi 11 ' 51 

The two catastrophes shown in Fig. 2] can affect the 
inherent structure network as shown in Fig. [5] We clas- 
sify catastrophes as the addition of either an internal or 
external edge to the network i 52 i 53 ' 54 i 55 An internal edge 
involves the addition of an edge between two existing 
nodes. This can occur on increasing p via a fold catastro- 
phe, creating a transition state and an index two saddle 
from an inflection point, or a cusp catastrophe, creating 
one transition state and two index two saddles from an 
index two saddle. Both these processes occur along the 




FIG. 5: (Colour online) Changes the inherent structure net- 
work that can occur on increasing p. (a) Addition of an in- 
ternal edge (shown in bold) connecting two existing nodes, 
El and E2. (b) Addition of an external edge connecting a 
new node N to an existing node E. (c) Addition of an exter- 
nal edge via a cusp. The initial transition state connects two 
permutational isomers, E and E' , forming a self-connection. 
The new node is inserted into this edge, forming connections 
to both permutational isomers. However, because we ignore 
ME/SCs, the overall effect on the network topology is equiv- 
alent to that in (b). 



boundaries between two basins of attraction. The new 
transition state may connect two permutational isomers 
of the same minimum, in which case it would form an SC, 
or it may connect two minima that are already connected 
by a transition state, forming a ME, and therefore have 
no effect on the network topology. 

An external edge involves the addition of a new node 
and one new edge. This can either be from a fold catas- 
trophe, an inflection point becoming a minimum and a 
transition state, or a cusp catastrophe, a transition state 
becoming a minimum and two transition states. Sym- 
metry restrictions on the cusp catastrophe affect the net- 
work as shown in Fig. [5Jc). The initial transition state 
connects two permutational isomers, therefore it forms 
an SC in the network and is not included in our repre- 
sentation. A cusp inserts a minimum into this edge, that 
is connected to both permutational isomers. As we only 
count one permutational isomer in the inherent structure 
network, the effect on the network topology is that one 
external edge is added. 

When a new minimum appears on the landscape, an 
external edge is formed to one existing node by the new 
transition state that was involved in the catastrophe. 
However, the change to the network topology can be more 
complicated. The basin of attraction associated with a 
new minimum generated by a fold catastrophe has a fi- 
nite area, and as such can cause reaction pathways from 
existing transition states to be rewired to connect to it 
if the existing reaction pathway passed through the re- 
gion of configuration space now assigned to the basin of 
the new minimum (Fig. [5]). The new node can therefore 
have initial degree greater than one, although there is 
only one new transition state. These rewired transition 
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(a) 





FIG. 6: (Colour online) The two major ways that the addition 
of an external edge (shown in bold) connecting new node N 
to existing node E can lead to an edge being rewired, thus 
increasing the initial degree of N. In (a), the edge is rewired 
from E, thus reducing its degree. In (b), the rewired edge was 
originally an ME, thus the degree of E is not affected. 



states may be rewired from MEs and as such appear to 
be new edges in terms of the network topology. Because 
such new edges connect to the new node, they are clas- 
sified as external edges. 

Although these are the most likely types of catastro- 
phe, there may be other, more complicated catastro- 
phes. However, such catastrophes involve higher sym- 
metries and are therefore rare. Furthermore, we do not 
search specifically for folds and cusps, simply for station- 
ary points on the PEL. Starting from p = 6, for which 
it is likely that all of the stationary points are known, 
we decrease p (increasing the range of the potential) and 
record when stationary points disappear from the sur- 
face, and when transition states rewire to connect to dif- 
ferent minima. This information is then inverted to give 
the details of the network's growth. While it is possi- 
ble for a minimum or transition state to appear as the 
range is increased, this type of event is difficult to detect 
and is expected to be relatively rare, so is not included. 
Rewiring of transition states is only considered when a 
minimum is involved in a catastrophe, as this is the most 
likely occasion for it to occur. 



III. MEASURING PREFERENTIAL 
ATTACHMENT 

In preferential attachment models, when a new node 
appears in the network, it is assumed to connect to an ex- 
isting node with probability II, depending on its degree k 
as II cx k a . The exponent a describes the strength of the 
preferential attachment. In the original Barabasi- Albert 
models a must be one to give a scale-free network^ In 
more complex models, allowing for example new edges 



between a pair of existing nodeS ) 57 i 58 ' 59 i 60 ' 61 i 62 other val- 
ues of a can give rise to scale-free networks. 

The exponent a can potentially be determined from 
time-resolved data for the growth of the network. How- 
ever, extracting the exponent is not straightforward for 
two reasons. Firstly, network growth is typically stochas- 
tic. A given preferential attachment rule can result in 
many different networks, and as such the exponent can- 
not be determined exactly. Secondly, at each time step 
the probability of connecting to a node with degree k 
is n = k a / y\ kf . The normalising sum in the de- 
nominator is time-dependent, and depends on the un- 
known exponent. Approximations have been used pre- 
viously, such as comparing the network at two different 
stages of its growth and assuming the normalising sum 
remains constant over this time period ) 58 ' 59 i 61 or assum- 
ing that the normalisation sum is proportional to the size 
of the network^ By contrast, we use an iterative, self- 
consistent method that we recently introduced^ In this 
method, we estimate a in order to evaluate the sum, then 
use the results to obtain a better estimate for a. 

The number of edges gained by all nodes with degree 
k at time t, Afc, depends on the number of nodes with 
degree k at time t, n k {k, t), and is given by 



Ak{k,t) = rik(k,t)IL(k,t) 



n k {k,t)k a _ n k (k,t)f(k) 



Eih(t) a c t {t) ' 

(2) 

where we define f(k) = k a and c t (t) = ^ ki(t) a . f(k) is 
obtained by summing both sides over all time steps, thus 
avoiding problems from steps where rik(k, t) — 0, giving 



f(k) 



S t Afc(fc,t) 

T,t n k{k,t)/Ct(t)' 



(3) 



a is then obtained from a linear least squares fit to 
log f(k) = alogfc, weighted by the number of terms 
in the sum in the denominator of Eq. [3] Logarith- 
mic binning over k reduces errors from degrees where 
53 t Afc(fc,t) — 0, which cannot be included in a log-log 
plot. It has been shown that the resulting exponent is 
independent of the initial estimate for a, and the method 
generally converges^ In tests on model networks grown 
with a known preferential attachment rule, the exponent 
obtained had an absolute error less than approximately 
0.13 for 0.5 < a < 1.5, with a much lower error for a 
closer to 1. 

We investigate the addition of internal and external 
edges independently. For internal edges, the probability 
of a pair of nodes gaining an edge is assumed to depend 
on the product of their degrees as H(kikj) cx (kikj) a . If 
multiple edges are not allowed in a network, when deter- 
mining a only those pairs of nodes that are not connected 
should be used in the normalization sum Ct(t). This 
modification has been seen to greatly improve the ex- 
ponent obtained, reducing the error from approximately 
20% to 1% in the case that a = lJ£ This approach can 
straightforwardly be applied to the evolution of the in- 
herent structure networks, the only unusual feature being 
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FIG. 7: (Colour online) Growth with p of the number of 
nodes, N, the number of edges, M, the average degree, 
(k) = 2M/N, and the initial degree of a new node, knew 
(binned over p) . 



that the networks grow as a function of the potential pa- 
rameter p, not time. 



IV. RESULTS 
A. Network properties 

The number of nodes and edges increase as the range 
parameter is increased towards p = 6, at which the net- 
work contains 1509 nodes and 20 691 edges. As shown 
in Fig. [3 the number of nodes increases roughly expo- 
nentially with p, i.e. the network grows faster as p ap- 
proaches 6. The number of edges also increases roughly 
exponentially with p. As cusp catastrophes occur at ei- 
ther transition states or index 2 saddles, and fold catas- 
trophes occur at inflection points, both of these will be 
more likely when the landscape is rougher (i.e. the net- 
work is bigger). In particular if dN/dp oc N exponential 
growth results. 

The number of edges increases faster than the number 
of nodes, so the average degree increases, i.e. the network 
is accelerating and becoming better connected. This is 
due to two factors. Firstly, the addition of internal edges. 
Secondly, the initial degree of a new node increases with 
p. If a node is added at higher p, it is likely to link 
to more existing nodes because there are more nodes in 
the network, and more edges to rewire. Some of these 
edges are rewired from MEs, thus increasing the average 
degree. 

The degree distribution has a similar form through- 
out the growth process, as shown in Fig.[5£a), indicating 
that the growth is fairly homogeneous in terms of network 
topology, and the network reaches a kind of steady state. 
Furthermore, the distribution of energies of the minima is 
approximately Gaussian throughout the growth process 
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FIG. 8: (Colour online) (a) Cumulative degree distribution 
(i.e. the number of nodes with degree greater than k) and (b) 
distribution of the energies of the minima at various values 
of p during the growth of the network. In (b), Gaussian best 
fits are also shown. In calculating these best fits, the global 
minimum is excluded as it has especially low energy, due to 
the particularly stable icosahedral structure of the cluster. 



(Fig. Efb)). Although the initial degree of 'newer' nodes 
(those added at higher p) is higher, the final degree at 
p = 6 is lower than for the older nodes (Fig. [Ha)). This 
"age" dependence of the degree is indicative of preferen- 
tial attachment and has been seen previously, for exam- 
ple in a protein interaction network.— It should be noted 
that the network is only large enough to study degree 
distributions for p > 3. In this region, the growth of the 
number of stationary points, and particularly the average 
degree, is also fairly steady (Fig. [7]). 



B. Preferential attachment 

The exponents that we obtained by applying the al- 
gorithm in Section Mil are shown in Table [U They are 
close to one, indicating that preferential attachment is 
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FIG. 9: (Colour online) Properties of a node at p = 6 as a 
function of the value of p when it was created, (a) Node de- 
gree, k, and (b) basin area, A. Individual points correspond 
to individual nodes and the solid line is a binned geometric 
average. Note that due to the exponential growth of the net- 
work, only 17 nodes, 1.1% of the total, appear at p < 3, where 
the correlation of the properties with p is weaker. 



present in the growth of the network. New edges are 
more likely to link to existing nodes with high degree, 
leading to a scale-free network. The exponent is slightly 
greater than one for external edges both for those asso- 
ciated with new transition states, and those associated 
with rewired MEs. When a new transition state appears 
on the PEL that gives rise to an external edge, the de- 
gree of the existing node to which it connects increases by 
one. Simultaneously, some edges are rewired from exist- 
ing nodes, increasing the initial degree of the new node. 
However, due to the spatial nature of the PEL, these 
rewired edges are likely to be rewired from the existing 
node that is connected to the new transition state (see 
Fig. [5]). In fact, on average, when such an external edge 
is added, the existing node at the end of the edge loses six 
edges to rewiring, whilst only gaining one. This opposes 
preferential attachment by taking edges from high-degree 
nodes. However, this effect is cancelled out by the pref- 
erential attachment associated with external edges that 



TABLE I: Exponents found for preferential attachment of 
external edges (II(fe) oc k a ) and internal edges (Jl{kikj) oc 
(kikj) a ). The number of data points used to obtain the expo- 
nents are also shown. External edges are further broken down 
into those that involve addition of a new transition state and 
those that involve rewiring a transition state. Rewired tran- 
sition states are only considered if they involve an increase in 
degree of the existing node, i.e. before the catastrophe, the 
transition state formed an ME. 

Exponent a Number 
External edges (new TS) 1.20 1 199 

External edges (rewired TS) 1.24 7 275 

Internal edges 0.99 12 274 



simultaneously appear due to the rewiring of MEs. Fur- 
thermore, the frequency with which external edges are 
added is lower than that for internal edges (Table [TJ) and 
the appearance of new minima increases the available 
basin boundary along which new internal edges can be 
generated. 

Internal edges do increase the degree of high-degree 
nodes, as expected for preferential attachment, and the 
exponent for this process is very close to one. Since these 
"time steps" dominate the network growth, the preferen- 
tial attachment associated with the addition of internal 
edges seems to be the main driving force for the scale- 
free topology of the final network. As internal edges 
are generated by catastrophes occurring along the basin 
boundaries, the larger, older minima are likely to garner 
more internal edges, hence the preferential attachment. 
Indeed, a correlation between degree and basin area is 
shown in Fig. [TOT aL 

To probe the connection to the basin area distribution 
in more detail, we have measured the areas of the basins 
of attraction directly. This task is achieved simply by 
choosing random points in configuration space and then 
performing a local minimization. The area of a basin is 
then proportional to the number of points that lead to 
the corresponding minimum after minimization. Newer 
minima have smaller basins of attraction in the final PEL 
at p = 6 (Fig. [Hb)). When new minima appear on the 
landscape they are likely to do so near to the edge of 
a basin of attraction. Near to the centre of the basin, 
the basin behaves harmonically, as quadratic terms dom- 
inate the Taylor expansion of the energy about the min- 
imum. Catastrophes are driven by higher-order anhar- 
monic terms, and so will only occur at sufficient distance 
from the minimum. 

At the exact moment at which the catastrophe occurs 
and the minimum appears in the network, the minimum 
and the new transition state also involved in the catas- 
trophe are at the same point. For the fold catastrophe, 
the basin associated with this minimum initially has zero 
extent in the direction of the transition state, but both in 
the opposite direction and all those directions perpendic- 
ular to this reaction coordinate, it is of finite extent. As 
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FIG. 10: (Colour online) (a) Degree dependence ol the area 
of the basin of attraction of a minimum for LJ13. Individual 
points correspond to individual minima and the solid line is 
a geometric average, (b) The cumulative basin area distribu- 
tion. In (b) a power-law form corresponding to p(A) ~ A~ 2 ' 2 
(the exponent of the cumulative distribution is greater by one) 
has also been plotted. 



at or near the boundaries of larger, existing basins. For 
this reason, the distribution of minima in configuration 
space will be inhomogeneous, with minima concentrated 
around the edges of the larger basins in a manner that 
is somewhat similar to the spatial distribution of hyper- 
spheres in the Apollonian packing (Fig. [2] 

Both the number of minima at a given p and the area 
of a basin (at p — 6) that appeared at p show an approxi- 
mately exponential dependence on p, i.e. N{p) ~ exp(i/p) 
and A(p) ~ exp(— ap). Hence we can derive an expres- 
sion for the area distribution. The number of basins with 
area A at p — 6 is given by 



n(A)dA 
n(A) 



dN , 
— dp 
dp 

i^exp(^p) 



dp 



dA 



v exp(t'p) 



Substituting for p via A(p) 



aexp(— ap) 
exp(— ap), we arrive at 



(4) 



n(A) ~ -A- X - V l a . 
a 



(5) 



This mechanism for generating a power law has been 
termed combination of exponentials by Newman^ The 
observed exponents are v w 1.1, a ~ 1.2, giving P(A) ~ 
A~ 19 . This result is in line with the observed approx- 
imate power-law tail to the basin area distribution for 
LJr3 (Fig. bVOTb)). The exponent is also close to that for 
the area distributions observed for liquids and for Apol- 
lonian packings. 12 For the latter the exponent is exactly 
—2 in the limit of high dimensionality. Interestingly, the 
analogies to the Apollonian packing go farther, as the 
power-law in the Apollonian case can also be considered 
to arise from a combination of exponentials. Namely, 
for the two-dimensional Apolloinian packing, it has been 
shown that both the number of disks and their area de- 
pend exponentially on the number of iteration steps used 
to generate the packing, and that, furthermore, new disks 
preferentially attach to those with high degreed 



p increases, the distance between the transition state and 
minimum will increase, albeit relatively slowly — close 
to the catastrophe this distance can be approximated by 
y/6AV/X, where AV is the bar rier, and A is the curvature 
at the minimum and transition state along the reaction 
path. 49 Therefore, the area of the new basin will proba- 
bly initially increase, but it is always likely to be smaller 
than the basins of the older minima between which it is 
sandwiched. Hence, as p increases and there are more 
minima on the landscape, the new minima will have in- 
creasingly small basins, as they are more constrained by 
the more closely spaced basins. 

Combining this effect with the growth of the network 
leads to a power-law area distribution. New basins that 
appear at higher p have smaller basin areas. Due to the 
exponential increase of N with p, there are many of these 
new, small basins. Furthermore, the nature of the catas- 
trophes means that the new basins are likely to appear 



V. CONCLUSION 

We have used techniques from network theory to study 
the evolution of a PEL on going from a simple form with 
a single minimum to a complex landscape with a large 
number of minima as the range of the underlying inter- 
atomic potential is varied. In particular, we have exam- 
ined how the inherent structure network, which provides 
a picture of the connectivity of the PEL, changes. New 
edges in this network are seen to preferentially attach 
to minima with a larger number of connections, helping 
to explain the scale-free nature of the inherent structure 
network. This preferential attachment is very close to 
linear, and can be understood by considering the under- 
lying catastrophes that lead to the appearance of new 
minima and transition states on the PEL. Each catas- 
trophe that occurs, as well as increasing the degree of 
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a minimum, also makes the surface rougher, leading to 
further catastrophes. These catastrophes are also more 
likely to occur at or near boundaries, therefore they are 
more likely to occur near large basins, which have high 
degree. 

Consideration of this growth process leads to further 
insight into the nature of the division of the PEL into 
basins of attraction. Previous work has shown that, sim- 
ilar to the Apollonian packing, PELs have power-law dis- 
tributions for the hyperareas of these basins of attraction, 
suggesting that the basins tile the landscape in a fractal- 



like manner. The current results are able to shed light on 
this puzzling character of the PEL. During the evolution 
of the PEL small new basins appear near the boundaries 
of larger existing boundaries, as in the generation of the 
Apollonian packing. The network grows exponentially, 
because catastrophes make the network larger and the 
surface rougher, enabling further catastrophes. The com- 
bination of this exponential growth and the smaller size 
of newer basins leads to a power-law area distribution, in 
agreement with observed results for PELs. 
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